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Abstract 


Two  canopy  windflow  models  (basic  model  and  alternative  model)  are  described  and  as¬ 
sessed  against  measurements  of  mean  wind  and  turbulence  in  a  wide  range  of  canopies 
(both  vegetative  and  urban).  The  windflow  models  are  based  on  a  first-order  turbulence 
closure  for  the  description  of  both  the  mean  wind  (IJ)  and  kinetic  energy  of  turbulence  (fc), 
viz  Reynolds5  equations  are  closed  using  eddy-viscosity  K  ~  A  A:1/2,  where  A  is  an  alge¬ 
braically  specified  turbulence  length  scale.  The  windflow  models  are  compared  with  eleven 
different  canopy  flows  (aeroelastic  plant  canopy,  corn  canopy;  square  and  staggered  arrays 
of  cubical  obstacles  at  three  different  frontal  area  densities  [0.0625,  0.16,  and  0  44],  square 
and  staggered  arrays  of  billboard-shaped  obstacles  at  a  frontal  area  density  of  0.16;  and, 
the  Tombstone  canopy).  Comparisons  of  model  predictions  with  these  different  canopies 
show  that  the  alternative  windflow  model  gives  the  best  overall  conformance  with  respect 
to  both  mean  flow  and  turbulence  profiles.  The  basic  windflow  model  was  found  to  under¬ 
predict  the  turbulence  kinetic  energy  within  the  canopy,  although  this  underprediction  does 
not  appear  to  affect  its  reasonably  good  predictions  of  the  mean  wind  speed  and  the  shear 
stress.  However,  both  the  basic  and  alternative  windflow  models  underpredict  severely  the 
magnitude  of  the  mean  wind  shear  at  the  canopy  top  for  for  “urban”  canopies  at  the  highest 
frontal  area  density  [0.44]  investigated.  It  is  suggested  that  the  alternative  windflow  model 
can  be  used  probably  unmodified  (viz.,  without  resorting  to  tuning  of  any  model  constants) 
for  the  simulation  of  a  wide  range  of  canopy  flows  provided  the  frontal  area  densities  of  the 
tkurban"  canopies  are  not  so  large  as  to  provoke  skimming  flow. 
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Executive  summary 

Introduction:  It  is  anticipated  that  Canadian  Forces  (CF)  in  the  foreseeable  future  will 
have  to  fight  in  or  protect  urban  areas,  whether  in  battle,  peace-making,  peace-keeping, 
or  counter-terrorist  operations.  The  increased  awareness  and  importance  accorded  by  the 
public  worldwide  and  their  governments  to  maintain  appropriate  defences  against  chemical 
and  biological  (CB)  agents  in  an  urban  (built-up)  environment,  the  prediction  of  casualties 
and  human  performance  degradation  resulting  from  such  releases,  and  the  development  of 
operational  procedures  and  regulations  to  control,  mitigate,  and  monitor  the  fate  of  CB 
agents  in  urban  areas  with  high  population  densities,  will  require  mathematical  modeling 
of  urban  wind  flows  and  dispersion.  In  this  regard,  it  should  be  noted  that  the  prediction  of 
urban  canopy  flow  is  in  principle  pre-requisite  to  or  co-requisite  with  the  prediction  of  the 
simpler  (but  itself  complex)  problem  of  scalar  (e.g  ,  CB  agent)  dispersion  within  a  canopy 

Results:  Two  canopy  wmdflow  models  (basic  model  and  alternative  model)  are  described 
and  assessed  against  measurements  of  mean  wind  and  turbulence  in  a  wide  range  of  canopies 
(both  vegetative  and  urban).  The  windflow  models  are  based  on  a  first-order  turbulence 
closure  for  the  description  of  both  the  mean  wind  ( U )  and  kinetic  energy  of  turbulence  ( k ), 
viz.  Reynolds’  equations  are  closed  using  eddy-viscosity  K  ~  Xkl// 2.  where  A  is  an  alge¬ 
braically  specified  turbulence  length  scale.  The  windflow  models  are  compared  with  eleven 
different  canopy  flows  (aeroelastic  plant  canopy;  corn  canopy;  square  and  staggered  arrays 
of  cubical  obstacles  at  three  different  frontal  area  densities  [0.0625,  0  16,  and  0.44],  square 
and  staggered  arrays  of  billboard-shaped  obstacles  at  a  frontal  area  density  of  0  16;  and, 
the  Tombstone  canopy).  Comparisons  of  model  predictions  with  these  different  canopies 
show  that  the  alternative  windflow  model  gives  the  best  overall  conformance  with  respect 
to  both  mean  flow  and  turbulence  profiles.  The  basic  windflow  model  was  found  to  under- 
predict  the  turbulence  kinetic  energy  within  the  canopy,  although  this  underprediction  does 
not  appear  to  affect  its  reasonably  good  predictions  of  the  mean  wind  speed  and  the  shear 
stress.  However,  both  the  basic  and  alternative  windflow  models  underpredict  severely  the 
magnitude  of  the  mean  wind  shear  at  the  canopy  top  for  for  “urban’’  canopies  at  the  highest 
frontal  area  density  [0.44]  investigated.  It  is  suggested  that  the  alternative  windflow  model 
can  be  used  probably  unmodified  (viz.,  without  resorting  to  tuning  of  any  model  constants) 
for  the  simulation  of  a  wide  range  of  canopy  flows  provided  the  frontal  area  densities  of  the 
“urban”  canopies  are  not  so  large  as  to  provoke  skimming  flow 

Significance  and  Future  Plans:  The  canopy  windflow  models  considered  in  this  report 
provide  all  the  statistics  of  the  wind  field  (i.e.,  spatial  variation  of  the  mean  wind  velocity, 
Reynolds  stresses,  and  viscous  dissipation  rate  over  the  urban  landscape)  required  as  the 
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input  to  drive  a  physically- based  model  for  the  prediction  of  urban  dispersion.  An  example 
of  such  a  physically- based  model  is  the  Lagrangian  Stochastic  (LS)  model.  This  model  re¬ 
quires  that  the  turbulent  velocity  field,  viz.  mean  velocity,  velocity  variances,  shear  stresses, 
and  the  Lagrangian  time  scale  (or,  equivalently,  the  viscous  dissipation  rate),  be  supplied 
as  the  external  inputs.  The  implication  of  the  differences  between  various  windflow  models 
for  mean  velocity,  velocity  variances,  and  viscous  dissipation  rate  should  be  tested  m  a 
Lagrangian  dispersion  model  to  determine  the  sensitivity  of  predicted  tracer  concentrations 
to  the  windflow  model  predictions  of  mean  flow  and  turbulence. 


Yee,  E.  and  Wilson.  J  D.,  2001  A  Comparison  Between  a  First-Order  Closure  Model  for 
Disturbed  Urban  Canopy  Flows  with  Observations  for  Mean  Wind  and  Turbulence.  DRES 
TR  2001-058,  Defence  Research  Establishment  Suffield. 
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Introduction 


It  is  anticipated  that  there  will  be  an  increasing  likelihood  that  Canadian  forces  in  the  fore¬ 
seeable  future  will  have  to  fight  in  or  protect  urban  areas,  whether  in  battle,  peace-making, 
peace-keeping,  or  counter-terrorist  operations.  The  increased  awareness  and  importance 
accorded  by  the  public  worldwide  and  their  governments  to  maintain  appropriate  defences 
against  chemical  and  biological  (CB)  agents  released  m  an  urban  environment,  the  predic¬ 
tion  of  casualties  and  human  performance  degradation  resulting  from  such  icleases,  and  the 
development  of  operational  procedures  and  regulations  to  control,  mitigate,  and  monitor 
the  fate  of  CB  agents  in  urban  areas  with  high  population  densities  will  require  the  de¬ 
velopment  of  physically-sound  and  fully-validated  models  for  the  disturbed  wind  flow  and 
dispersion  within  an  arbitrary  domain  inside  a  city  (or,  other  built-up  area). 

In  principle,  particle  trajectory  (Lagrangian)  dispersion  models  can  be  used  to  calculate 
dispersion  in  an  urban  environment  provided  the  spatial  pattern  of  mean  wind  and  tur¬ 
bulence  statistics  over  the  urban  landscape  is  known.  Consequently,  an  understanding  of 
the  complex  flow  patterns  which  result  from  the  interaction  of  the  wind  with  clusters  of 
buildings  is  required  in  order  to  predict  the  dispersion  of  harmful  materials  released  in  an 
urban  environment  To  assist  in  this  understanding,  there  have  been  attempts  to  apply 
canopy  flow  models  developed  originally  for  vegetation  (e.g.,  forests,  crops,  etc.)  to  the 
urban  landscape.  This  is  a  natural  development  since  vegetative  canopy  flow  models,  capa¬ 
ble  of  reproducing  many  features  of  the  velocity  statistics  needed  in  Lagrangian  stochastic 
(LS)  models,  have  been  developed  over  the  past  thirty  years  with  varying  complexity,  lim¬ 
itations,  and  successes.  These  models  range  from  simple  algebraic  mixing  length  models 
[1] ,  [2]  to  one-  and  two-equation  prognostic  models  for  the  determination  of  the  eddy  viscos¬ 
ity  [3], [4]  Second-order  local  closure  models  [5]-[7],  and  still  higher-order  models  [8]  have 
been  developed  for  plant  canopy  flows  Although  limitations  of  first-order  closure  models 
for  canopy  flows  are  well- recognized  [9],  they  nevertheless  are  still  widely  employed.  This  is 
because  these  models  provide  both  the  simplicity  and  flexibility  for  describing  many  of  the 
relevant  turbulent  processes  in  canopy  flows,  while  enabling  the  solution  of  the  resulting 
model  equations  to  be  computationally  tractable  (viz.,  second-  and  higher-order  models  are 
much  more  computationally  demanding  than  first-order  closure) 

The  primary  objective  of  this  study  is  to  describe  and  investigate  whether  a  simple  first- 
order  closure  model  for  the  prediction  of  the  flow  within  and  above  a  canopy  (either  vege¬ 
tative  or  urban)  yields  mean  velocity  and  second-order  velocity  statistics  that  are  in  close 
enough  agreement  with  measurements  so  as  to  provide  adequate  input  information  for  use 
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m  Lagrangian  dispersion  models.  For  this  purpose,  a  number  of  canopy  experiments  (some 
vegetative,  but  mostly  urban)  have  been  identified  for  use  in  model  validation.  We  restrict 
our  attention  to  experiments  providing  measurements  of  both  the  mean  flow  and  the  tur¬ 
bulence  (turbulence  kinetic  energy,  shear  stress,  and  normal  stresses),  and  depend  heavily 
on  high-quality  data  from  wind  tunnel  and  water  channel  studies 
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Windflow  model 


A  detailed  explanation  of  the  model  equations  and  numerical  procedure  is  given  in  Yee  and 
Wilson  [10],  and  so  we  shall  give  only  a  brief  description  here.  In  the  so-called  micro-scale 
(horizontal  scales  less  than  2  km)  resolution  of  atmospheric  flow  through  an  urban  canopy, 
we  solve  the  continuity  equation  (mass  conservation)  and  the  Reynolds-averaged  Navier 
Stokes  equations  (momentum  conservation)  for  steady  and  incompressible  fluid  flow  under 
adiabatic  conditions  (neutrally-stratified  atmosphere)* 

Continuity  equation . 


ou  ov  d\v 

dx  +  dy  dz 


(1) 


Reynolds- aver  aged  Navier -Stokes  equations : 


0  —77 


0 


d 


—  (U2  +  a'2)  -f  —  (UV  4-  uV)  +  —  (UW  +  u'w ')  =  f(V  -  VG)  - 


dx 


dy 


dx 


-  CdAU \flP-  +  V?  +  W2. 


(2) 


C\  r\  n  rj  rj 

—  (UV  +  ^P)  +  —(V2  +  V' 2)  +  —  (VW  +  i/w7)  =  f(UG  -U)-  — 

-  CdAV \JlP  +  V 2  +  W2,  (3) 

and 

O  Q  rj  T-) 

—  {UW  +  W)  +  —  {VW  +  Vn?)  +  —  (W2  +  W*)  =  -  — 
dx  dy  dz  dz 

-  CdAW  \J U2  +  V2  +  l^2.  (4) 

corresponding  to  the  streamwise,  cross-stream,  and  vertical  mean  momentum  balance,  re¬ 
spectively  Here,  x  =  (xi,  X2,  X3)  =  (x,y,z)  is  the  position  vector  with  the  x-coordmate  in 
the  mean  streamwise  direction,  the  ^-coordinate  in  the  horizontal  cross-stream  direction, 
and  the  2-coordmate  normal  to  the  ground,  U  =  (U 1,  U2 ,  C/3)  =  ( U ,  V,  W)  is  the  mean  veloc¬ 
ity  vector  with  components  in  the  x-,  7/-,  and  ^-coordinates,  respectively;  u'  =  (ix'l5  u'2,  u'3)  = 
(u\v\wr)  is  the  velocity  fluctuation  vector  (departures  from  the  mean  velocity  vector)  in 
the  x-,  y-7  and  c-directions,  respectively,  P  is  the  mean  kinematic  pressure;  /  is  the  Coriolis 
parameter  defined  as  /  =  2usin<£,  where  u ;  is  the  magnitude  of  the  angular  velocity  of  the 
Earth’s  rotation,  and  d>  is  the  latitude;  C d  is  the  bulk  drag  coefficient  assigned  to  either  a 
single  canopy  element  (obstacle)  or  a  whole  group  of  canopy  elements;  A  —  A(x,y,z)  is  the 
frontal  area  per  unit  volume  (m2  m~3)  of  space  occupied  by  either  a  single  canopy  element 
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or  a  whole  group  of  canopy  elements;  and,  an  over  bar  is  used  to  represent  a  Reynolds  av¬ 
erage.  In  Equations  (2)  and  (3),  the  large-scale  background  horizontal  kinematic  pressure 
gradient  in  the  Planetary  Boundary  Layer  (PBL)  has  been  parameterized  in  terms  of  a 
nominal  geostrophic  wind  aloft,  so  fUc  —  and  f^G  —  ^ r ,  where  ( Ug,Vg )  are  the 

components  of  the  geostrophic  wind  aloft. 


The  Reynolds-averaged  Navier-Stokes  equations  (2)  to  (4)  are  not  closed  unless  a  model  is 
provided  which  ties  the  kinematic  (specific)  Reynolds  stress  tensor  u[u '  (i,j  =  1,2,  3)  to  the 
global  history  of  the  mean  fluid  velocity  Uz  in  a  physically  consistent  fashion.  Herein  lies 
the  fundamental  problem  of  turbulence.  In  order  to  compute  all  the  mean-flow  properties 
of  the  turbulent  flow  under  consideration,  we  need  a  prescription  for  computing  u'u'  The 
specific  Reynolds  stresses  appear  in  Equations  (2)  to  (4)  as  spatial  gradients  and,  physically, 
can  be  interpreted  as  forces  "felt  by 7  the  mean  flow.  The  simplest  method  for  relating  the 
Reynolds  stresses  to  the  mean  velocity  is  to  use  the  Boussinesq  linear  and  isotropic  eddy 
viscosity  approximation 


— u'u'  —  2  KSn  —  -kSt 

l  J  LJ  o  L. 


(5) 


where  K  is  the  eddy  viscosity,  k  is  the  turbulence  kinetic  energy  (TKE),  Sl3  =  \{j^r  + 
)  is  the  mean  strain-rate  tensor,  and  StJ  is  the  Kronecker  delta  function.  As  a  further 
simplification,  we  note  that  in  most  disturbed  thin  shear-layer  flows  (i.e  ,  flows  wherein  shear 
strains  normal  strains),  the  normal  stress  gradients  play  a  much  smaller  role  than  the 
shear  stress  gradients.  In  view  of  this,  we  will  solve  simplified  £/,  V,  lE-momentum  equations, 
that  represent  only  what  are  (according  to  experience)  the  dominant  terms  (advection  by 
the  mean  flow,  pressure  gradient,  drag  on  obstacles,  and  the  divergence  of  the  shear  stresses 
[turbulent  momentum  fluxes]  in  the  three  coordinate  directions).  With  these  simplifications, 
the  mean  momentum  equations  can  be  written  explicitly  as 


d 

dx 
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=  f(uG  cdAVy/u*Tv*Tm,  (7) 
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The  inclusion  of  the  normal  stress  gradients  m  Equations  (6)  to  (8)  through  the  Boussmesq 
eddy- viscosity  approximation  of  Equation  (5)  does  have  the  purely  computational  advantage 
of  providing  ‘stabilizing"  diffusion  terms.  For  simplicity,  we  have  neglected  the  normal  stress 
gradients,  but  to  ensure  numerical  stability,  we  added  diffusion  terms  m  Equations  (6)  to  (8) 
where  diffusion  was  otherwise  absent  (through  the  neglect  of  the  normal  stress  gradients). 
We  regard  these  as  purely  artificial,  and  signify  the  artificial  diffusivity  as  Ka,  which  is  very 
small.  Hence,  Ka  is  used  for  no  other  purpose  than  to  stabilise  the  numerical  solution  of 
Equations  (6)  to  (8). 


The  only  unspecified  quantity  in  Equations  (6)  to  (8)  is  the  (true)  eddy  viscosity.  K  The 
latter  may  always  be  regarded  as  the  product  ( K  —  Ar)  of  a  turbulence  length  scale  (A)  and 
velocity  scale  (r).  Its  specification,  propitious  or  otherwise,  is  often  the  key  to  the  success 
of  a  flow  simulation  In  our  case,  the  eddy  viscosity  is  estimated  as 

K  —  A(x,  y,  z)  yj cek(x .  y,  z)  (9) 


Here,  the  turbulence  length  scale  A  is  a  function  of  x,  y,  and  z  for  an  inhomogenous  urban 
canopy  In  this  specification  for  K,  A;1/2  has  been  used  as  the  turbulence  velocity  scale.  The 
constant,  ce,  in  Equation  (9)  is  the  equilibrium  shear  stress/TKE  ratio  (u20/A:o)  immedi¬ 
ately  above  the  canopy  It  should  be  noted  that  this  constant,  when  evaluated  in  the  local 
equilibrium  region  of  a  free-shear  or  boundary-layer  flow,  varies  little  from  one  flow  config¬ 
uration  to  another  In  particular,  various  measurements  which  have  been  summarized  by 
Townsend  [11]  indicate  that  for  boundary  layers,  wakes,  and  mixing  layers,  ce  is  nearly  the 
same  and  varies  between  about  0.2  to  0.3  The  TKE,  A:(x,y,  z),  is  calculated  as  a  dynamic 
variable  from  the  following  approximate  transport  equation. 


d_ 

dx 


=  K 


d  (  dk 


+ h  (Wk  - 


dU  dV 
dy  dx 
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ou_  owv 

dz  +  dx  ) 
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dv  dwy 
+  ~dy  J 


dz 


e. 


(10) 


where  P w  is  the  wake  production  term  and  e  is  the  viscous  dissipation  of  TKE.  The  constant 
li  in  Equation  (10)  represents  the  ratio  of  the  effective  eddy  diffusivity  for  TKE  to  the  eddy 
viscosity  (i.e..  inverse  of  a  turbulent  Schmidt  number).  For  simplicity,  in  the  basic  version 
of  this  wind  flow  model,  it  is  assumed  that  Pw  =  0  From  this  simplified  perspective,  we 
regard  k  as  excluding  the  TKE  residing  in  the  wake-scales,  in  which  case  the  gain  to  k 
from  the  mean  kinetic  energy  (MKE)  is  excluded.  The  TKE  dissipation  (i.e.,  sink  — e)  in 
Equation  (10)  is  then  specified  as  e  —  max(ecc,  efd)  where 


(cek)3/2 


and  €fd  =  a,CdAkyftj*  +  V2  +  W2. 


(11) 
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The  term  ecc  is  the  standard  parameterization  for  viscous  dissipation,  the  only  TKE  sink 
above  the  canopy  (i.e.,  €cc  balances  shear  production  of  TKE,  in  the  local  equilibrium  layer 
well  above  the  canopy).  Our  specification  for  (“form  drag"),  which  converts  resolved 
TKE  to  the  small,  rapidly  dissipated  “wake-scales",  can  be  justified  on  physical  grounds. 
It  is  expected  that  a*  is  an  order  one  closure  constant,  but  its  exact  value  depends  on  the 
partitioning  of  k  between  the  three  velocity  components  (u' ,  v\  wf ).  Clearly,  the  particular 
canopy  geometry  plays  a  role  in  the  energy  partition. 


It  remains  to  specify  the  turbulence  length  scale  A.  Ideally,  this  turbulence  length  scale 
should  be  determined  as  a  dynamic  variable  for  the  case  of  inhomogeneous  canopy  flows 
but,  for  simplicity,  we  choose  instead  to  specify  this  quantity  algebraically.  In  this  sense, 
the  current  windfiow  model  is  incomplete  since  the  turbulence  length  scale  is  not  provided 
automatically  by  the  model  proper,  but  needs  to  be  specified  in  terms  of  some  relevant  flow 
characteristics  of  the  canopy  To  this  purpose,  the  turbulence  length  scale  is  modeled  as 
A  —  max(A2,  A0)  where 


_1_ 

K 


l  l 

- h  — 


and 


J_  -  1  1 

Ao  d)  ^OO 


(12) 


The  "outer"  scale  A0  recognizes  the  displacement  (of  the  equilibrium  structure)  due  to  the 
canopy  (through  the  zero  plane  displacement  d)  and,  through  the  imposed  (i  e.,  maximum 
permitted)  value  the  restricted  depth  of  the  PBL.  In  the  case  of  a  full  PBL  simula¬ 
tion,  Zoo  should  limit  the  length  scale  to  the  PBL  depth  which  is  often  estimated  [12]  as 
approximately 


_  0.0004C/g 

~  7 


(13) 


The  “inner"  scale  A*  recognizes  both  the  limitation  on  eddy  size  by  proximity  to  the  ground 
and  the  presence  of  the  canopy  The  influence  of  the  canopy  on  the  length  scale  is  embodied 
in  Ac,  the  canopy  shear  length  scale  which  is  taken  as 


Ac(x,y) 


ckl/2(x,y,z  =  hc ) 


dS 

dz 


z=hc 


(14) 


where  S  =  (U2  +  V 2)02  is  the  horizontal  windspeed  and  c  is  an  order  one  [viz.,  c  ~  0(1)} 
closure  constant.  The  canopy  shear  length  scale  has  been  shown  [13]  to  play  a  critical  role 
in  the  dynamical  process  that  produces  the  dominant  canopy  eddies. 


The  basic  windfiow  model  described  above  has  three  closure  constants  that  need  to  be 
specified,  namely,  the  inverse  of  the  turbulent  Schmidt  number,  /i  [Equation  (10)],  the 
coefficient,  a*,  specifying  the  form  drag  dissipation  [Equation  (11)];  and,  the  constant  of 
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proportionality,  c,  that  relates  the  canopy  shear  length  scale  to  turbulence  velocity  scale 
and  mean  shear  at  canopy  height  [Equation  (14)].  If  our  theory  were  exact,  we  could  set  the 
values  of  these  coefficients  from  first  principles  Unfortunately,  the  theory  is  not  exact,  but 
rather  a  model  developed  mainly  on  the  strength  of  dimensional  analysis  Consequently, 
the  best  we  could  do  is  to  set  the  values  of  the  closure  coefficients  to  assure  agreement 
with  observed  properties  of  canopy  flows.  To  this  purpose,  Wilson  et  al.  [14]  calibrated  the 
basic  wmdflow  model  against  windflow  observations  obtained  upwind  from  the  ridge  in  the 
Furry  lullflow  experiments  [15]  The  optimization  of  the  closure  constants  with  respect  to 
the  Furry  Hill  upwind  observations  yielded  the  following  values  for  the  closure  coefficients* 
c  =  1.0,  a*  =  1.0,  and  (i  —  02 

In  addition  to  the  basic  windflow  model  described  above,  we  will  consider  also  an  alter¬ 
native  variant  of  this  model.  Although  three  closure  constants  have  been  specified  in  the 
basic  wmdflow  model  by  calibrating  the  closure  scheme  against  the  Furry  Hill  data,  m  the 
alternative  variant  of  this  model  we  consider  the  removal  of  two  of  these  closure  constants. 
According  to  Launder  and  Spalding  [16,17],  the  specification  fi  =  1  0  (equality  of  eddy 
viscosity  and  eddy  diffusivity  for  TKE  or,  equivalently,  the  turbulent  Schmidt  number  for 
TKE  is  unity)  is  optimal  for  free  shear  flows  (plane  jets  and  mixing  layers);  and,  also  for 
boundary-layer  flows,  //  should  have  a  value  that  is  approximately  unity.  In  view  of  the 
similarity  of  canopy  flows  and  plane  mixing  layer  flows  [13],  it  is  surprising  that  ji  =  0.2  is 
optimal  based  on  comparison  with  the  Furry  Hill  upwind  observations  It  is  possible  that 
our  specification  fi  =  0.2  in  the  basic  model  is  due  to  deficiencies  in  our  representation  of 
other  terms  in  the  transport  equation  for  TKE.  particularly  in  the  wake  production  and 
form  drag  dissipation  terms.  In  view  of  this  possibility,  in  the  alternative  windflow  model, 
we  adopt  the  standard  recommendation  that  ft  =  1  [16,17]  and  consider  more  carefully  the 
parameterization  of  the  wake  production  and  form  drag  dissipation  terms. 

In  the  basic  windflow  model,  the  wake  production  term  Pw  was  assumed  to  be  zero  (viz., 
the  MKE  conversion  term  is  omitted)  with  the  effective  gam  to  the  TKE  from  the  MKE  and 
loss  from  the  TKE  to  smaller  dissipative  scales  being  parameterized  in  terms  of  a  viscous 
dissipation  of  the  form  e  —  max(6cc,  6fj),  with  the  closure  constant  a*  associated  with  our 
rather  inelegant  representation  of  the  form  drag  dissipation  [cf.  Equation  (11)].  In  view  of 
the  solid  bluff  obstacles  comprising  an  urban  canopy,  it  is  reasonable  for  these  obstacles 
(with  characteristic  dimensions  comparable  to  the  canopy  height)  to  transfer  fc,  not  to 
irrelevant  fine-scale  wake  scales  that  are  rapidly  dissipated  but,  rather,  to  larger  eddies  in 
the  energetic  sub-range  and/or  the  inertial-convective  sub-range  that  will  be  involved  in  a 
further  dyanamical  role  in  the  canopy  processes  To  parameterize  the  wake  production,  we 
make  the  assumption  that  the  partition  of  the  resulting  TKE  between  velocity  components 
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is  isotropic  at  the  wake  scales  (approximately  or  better)  and  write 

P™  =  CdA(U2  +  V2  +  W2)3/2^,  (15) 

as  the  wake  production  contribution  to  u\ufd  (creation  of  u'u'  as  the  mean  flow  does  work 
against  canopy  drag).  In  consequence,  the  wake  production  of  TKE  is 

Pw  =  \P™  =  \cdA(U 2  +  y2  +  IV'2)372.  (16) 

In  summary,  in  our  alternative  wmdfiow  model  we  retain  wake  production  in  the  form 
parameterized  by  Equation  (16)  as  a  gain  term  in  the  TKE  transport  equation  because 
urban  canopies  tend  to  be  dominated  by  large-scale  structures  (buildings)  with  dimensions 
and  separations  of  the  same  order  as  the  height  of  the  canopy. 


In  retaining  the  wake  production  term,  we  consider  also  a  more  careful  treatment  of  the 
process  of  dissipation  of  this  wake  TKE  in  our  alternative  windflow  model.  The  form  drag 
dissipation  accounts  for  the  conversion  of  resolvable  TKE  to  heat  in  the  obstacle  wakes  by 
working  of  the  turbulent  flow  against  obstacle  drag.  The  instantaneous  drag  force  acting 
in  the  ?-th  coordinate  direction  can  be  written  as  Ft  =  — C'dA||u||ul  ( i  =  1,2.3),  where 
u  =  (rti,  U2,  U's)  =  [U\  -b  u[ ,  U2  T  1*2,  C/3  +  U3)  is  the  total  instantaneous  velocity  vector  and 
||  *  ||  denotes  the  Euclidean  norm  of  a  vector  The  dissipation  of  the  specific  Reynolds  stress 
u'u'  due  to  form  drag  (which  is  denoted  by  e\Jd)  arises  from  the  working  of  the  turbulence 
against  the  fluctuating  drag  force  (F')  and  can  be  obtained  from  the  following  expression: 


4S  =  u[f;  +  u'jF; 


=  cdA(u,m<  +  tf.HUK  +  2||U||u>;). 


(17) 


Now.  note  that  the  Euclidean  norms  of  the  total  instantaneous  velocity  vector,  the  mean  ve¬ 
locity  vector,  and  the  fluctuating  velocity  vector  are  related  as  ||u||  =  ||U||  +  ||u,||.  Inserting 
this  relationship  in  Equation  (17)  and  neglecting  all  terms  involving  third-order  moments 
of  the  velocity  fluctuations,  we  obtain  the  following  approximate  expression  for  . 

eyd  «  CdA(u3WR  +  UtWWj  +  2||U||<^)  (18) 


As  it  currently  stands,  Equation  (18)  cannot  be  used  in  our  alternative  windflow  model 
because  we  need  some  method  to  estimate  the  correlation  ||u'||n'.  To  this  purpose,  we  first 
observe  that  ||u'||u'  can  be  expressed  as 


u  w  = 


null  1  + 


2  Uku'k 

l|U||2  '  ||U|| 


+ 


ukuk 


1/2 


I|U|| 


u,. 


(19) 
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Now  apply  a  binomial  series  expansion  to  the  square  root  term  m  Equation  (19)  and 
ensemble- average  the  resulting  expression,  neglecting  all  terms  involving  velocity  moments 
higher  than  second-order,  to  get 


Ukufku[ 

HUH 


(20) 


Finally,  inserting  Equation  (20)  into  Equation  (18)  gives  the  following  practical  estimate 

*3 . 

fd 

eyd «  CM  [2\\V\\M  +  U,Uk^r  +  uxuk- 


for  elJ  • 


«  ■  rrrr  \ 

nun +  1  nun ) 


(21) 


In  view  of  Equation  (21),  the  form  drag  dissipation  of  TKE  can  be  obtained  from 


£fd  —  Tr^fci  “  CdA  [  2 1 1 U j | A;  +  UxUk 


ukui 

nun 


(22) 


The  parameterization  for  era  in  Equation  (22)  [which  we  use  in  the  alternative  windflow 
model]  does  not  involve  any  closure  constants,  in  contrast  to  the  parameterization  of  Equa¬ 
tion  (11)  used  in  the  basic  windflow  model. 


In  summary,  the  alternative  windflow  model  eliminates  two  closure  constants  used  in  the 
basic  windflow  model  Firstly,  the  turbulent  Schmidt  number  fi~x  is  set  exactly  to  unity, 
which  is  consistent  with  others7  shear  layer  simulations  where  the  eddy  viscosity  is  assumed 
to  be  identical  to  the  turbulent  diffusivity  for  TKE  Secondly,  Equation  (22)  is  used  as  the 
parameterization  for  the  form  drag  dissipation  6fcj,  eliminating  the  closure  constant  a *  The 
remaining  closure  constant  c  (which  arises  as  the  proportionality  constant  in  specification 
of  the  canopy  shear  length  scale)  retains  the  same  value  that  it  has  m  the  basic  windflow 
model  (namely,  c  =  1). 
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Experimental  data  for  model  evaluation 


The  ideal  way  to  assess  the  accuracy  of  the  basic  and  alternative  canopy  windflow  models 
is  to  compare  their  performances  for  a  broad  range  of  canopy  test  flows.  We  restrict  our 
attention  to  experiments  providing  measurements  of  both  the  mean  flow  and  the  turbulence 
(shear  stress.  TKE.  and/or  normal  stresses),  and  depend  heavily  on  data  from  controlled 
experiments  in  a  wind  tunnel  or  water  channel.  A  number  of  vegetative  and  urban  canopies 
has  been  chosen  for  this  model  evaluation,  but  the  emphasis  has  been  placed  on  canopies 
composed  of  solid  bluff  obstacles  (ie.,  “urban"  canopies).  The  canopies  considered  here 
are  uniform  In  an  ideal  assessment  of  the  accuracy  of  a  model,  extraneous  errors  (e.g., 
measurement  error,  numerical  error,  discrepancies  in  the  boundary  conditions,  etc.)  need 
to  be  shown  to  be  small  so  that  the  discrepancy  between  measured  and  calculated  flow 
properties  arise  only  from  inaccuracies  of  the  model  (i.e.,  model  error),  the  test  flows  should 
be  different  than  those  used  in  the  development  of  the  model,  and,  the  tests  should  not  be 
performed  by  the  model's  developers  Even  though  the  model  evaluation  reported  here  does 
not  conform  to  this  ideal  m  all  respects,  it  is  nevertheless  hoped  that  the  results  are  useful  in 
helping  the  reader  assess  the  accuracy  and  utility  of  the  proposed  canopy  windflow  models. 


For  each  experiment,  we  determined  a  bulk  drag  coefficient,  Cd,  using  the  information 
from  the  measured  mean  wind  and  shear  stress  within  the  canopy  To  this  end,  Cd  was 
determined  from 


r(hc) 

?.2 


=  CdAhr 


(23) 


where  r(hc)  is  the  shear  stress  at  the  canopy  top  (z  =  hc)  [i.e..  r  =  —  u'w'}  and  u *0 
is  the  reference  friction  velocity  derived  from  the  shear  stress  at  canopy  top.  In  writing 
Equation  (23),  we  have  explicitly  assumed  that  the  canopy  is  uniform.  Furthermore,  in  our 
model  simulations  as  well  as  in  the  presentation  of  the  experimental  data.  u*0  is  used  as  the 
reference  velocity  scale  for  the  normalization  of  all  flow  quantities.  From  observations  of  au , 
<rv,  and  (Tw  (standard  deviations  of  velocity  fluctuations  in  the  streamwise.  cross-stream, 
and  vertical  directions,  respectively)  above  the  canopy,  we  determined  equilibrium  ratios 
Cu  =  cv  =.  cTyjxt^Q  and  cw  —  oru?/^*o  and 


Ce 


2 

cl  +  C2  +  c2w 


(24) 


For  the  canopy  flows  in  wind  tunnel  and  water  channel  experiments,  there  was  a  vertical 
gradient  in  the  shear  stress  above  the  canopy  (balanced  by  advection).  Normally,  in  these 
kinds  of  experiments  the  inertial  sublayer,  where  the  shear  stress  is  constant  and  where  the 
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mean  velocity  adheres  to  the  logarithmic  law-of- the- wall  profile,  is  “squeezed”  out  between 
the  roughness  sublayer  and  outer  layer  owing  to  the  large  value  of  hc/S  (where  5  is  the 
boundary-layer  depth)  compared  with  the  atmosphere  In  consequence,  m  order  to  properly 
account  for  stress  and  TKE  gradients  at  z  >  hc  in  the  canopy  flow  experiments  conducted  m 
a  wind  tunnel  or  water  channel,  we  incorporated  an  effective  background  kinematic  pressure 
gradient  which  was  assigned  from  the  measured  shear  stress  gradient  above  the  canopy,  viz. 


1  0Pe  _  d  f  r 

u*o  dx/hc  dz/hc\ulQ 


(25) 


where  P€  is  the  effective  kinematic  pressure.  The  value  for  the  effective  kinematic  pressure 
gradient,  determined  in  Equation  (25)  for  the  experiment,  was  imposed  in  the  simulations 
of  the  experiment  The  canopies  used  in  this  model  evaluation  are  described  briefly  below. 
The  first  two  canopies  are  vegetative  canopies  composed  of  flexible  elements;  the  remaining 
nine  canopies  correspond  to  urban  canopies  composed  of  solid  bluff  obstacles 


Brunet  et  al.  [18]  described  a  wind  tunnel  “wheat”  canopy  in  detail  that  consisted  of  an  array 
of  stalks  constructed  from  monofilament  nylon  fishing  line  of  diameter  0.25  mm,  spaced  on  a 
uniform  square  grid  of  side  5  mm.  Each  stalk  was  nominally  50  mm  long  but,  because  of  the 
natural  tendency  of  the  stalk  to  bend,  the  mean  canopy  height  hc  was  only  47  mm  In  the 
Furry  Hill  experiment  [19],  this  canopy  was  placed  upwind  and  over  a  Witch  of  Agnesi  ridge. 
The  parameters  of  this  canopy,  relevant  for  the  present  work,  are  summarized  in  Table  1. 
For  the  specification  of  av  (which  was  not  measured)  and  the  TKE,  the  approximation  that 
au  —  av  was  used. 


Table  1.  Parameters  of  the  Furry  Hill  experiment . 


Property 

Symbol 

Value 

Canopy  height 

hc 

0.047  m 

Bulk  drag  parameter 

CdAhc 

0.32 

Frontal  area  per  unit  volume 

A 

10.0  in-1 

Friction  velocity 

0.975  m  s"1 

Mean  wind  at  z  =  hc 

P  {hc)  I'm*  o 

3.7 

Equilibrium  standard  deviations 

(<7U,  O’ vi  Ow) / XL 

2.2,  2.2,  1.25 

TKE  at  z  =  hc 

k(hc)/ul  o 

5  62 

Effective  pressure  gradient 

(hc/ul0)dPe/dx 

-0.16 

Displacement  height 

d 

0.0333  m 

The  Elora  corn  canopy  was  described  by  Wilson  et  al.  [20]  who  presented  profiles  of  mean 
wind  and  turbulence  statistics  within  (but  not  above)  a  mature  corn  canopy  at  Elora, 
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Ontario.  The  measurements  of  the  flow  statistics  were  made  using  servo-controlled,  split- 
film,  heat-transfer  anemometers.  Some  relevant  parameters  for  this  canopy  are  presented  in 
Table  2.  A  deep  constant  stress  layer  was  used  as  the  background  flow  for  this  simulation 

Table  2.  Parameters  of  the  Elora  com  canopy  experiment. 


Property 

Symbol 

Value 

Canopy  height 

hc 

2.21  m 

Bulk  drag  parameter 

CdAhc 

0.79 

Mean  wind  at  2  =  hc 

r {hc) / u* 0 

3.04 

Equilibrium  standard  deviations 

(cTuj  &vi  O’u1)/ u*0 

2.06,  1  65.  1  13 

TKE  at  z  —  hc 

k(hc)/u2,  o 

4.12 

Displacement  height 

d 

15m 

Macdonald  et  al.  [21]  reported  profiles  of  mean  wind  and  turbulence  statistics  in  square 
arrays  of  cubical  obstacles  with  frontal  area  indices  Xp  (frontal  area  of  obstacle  normalized 
by  the  lot  area  within  which  a  single  obstacle  sits)  of  0.0625,  0  16,  and  0  44  corresponding, 
respectively,  to  isolated  roughness  flow,  wake  interference  flow,  and  skimming  flow.  In 
these  square  (m-line)  arrays,  there  were  unobstructed  streets  between  the  obstacles  in  the 
direction  of  the  approaching  wind  These  expeiiinents  were  conducted  m  a  water  channel 
using  a  SonTek  MicroADV  acoustic  Doppler  velocimeter  to  measure  the  flowr  velocities.  The 
parameters  for  these  square  arrays  are  given  in  Table  3  In  all  these  obstacles  arrays,  the 
frontal  area  density  A  can  be  determined  from  A  =  Xp/hc. 

Macdonald  et  al.  [21]  provided  data  on  mean  wind  and  turbulence  within  and  above  stag¬ 
gered  arrays  of  cubical  obstacles  with  frontal  area  indices  Xp  =  0.0625.  0  16,  and  0  44. 
In  these  uniform  obstacles  arrays,  alternate  rows  of  obstacles  (cubes)  were  shifted  by  half 
the  center-to-center  spacing  between  the  obstacles  For  the  staggered  arrays  of  cubes  with 
Xp  =  0.0625,  0.16,  and  0  44,  the  parameters  relevant  to  the  model  simulations  can  be 
summarized  in  Table  4 

The  last  set  of  obstacle  arrays  selected  are  square  and  staggered  arrays  of  billboard-shaped 
obstacles.  Unlike  cubical  obstacles,  the  billboard-shaped  obstacles  have  small  depth  with 
the  result  that  the  obstacles  are  effectively  two-dimensional  (approximately  or  better)  Mac¬ 
donald  et  al.  [21]  measured  the  mean  wind  and  turbulence  statistics  within  and  above  square 
and  staggered  arrays  of  billboard-shaped  obstacles  with  a  frontal  area  index  Xp  —  0.16.  The 
obstacles,  machined  from  aluminum  angle  sections,  were  50  mm  wide,  3  mm  thick,  and  50 
mm  high.  The  relevant  parameters  for  the  square  and  staggered  arrays  of  obstacles  are 
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given  in  Table  5. 

Table  3.  Parameters  of  the  square  arrays  of  cubical  obstacle  experiments. 


Property 

Symbol 

Value 

Frontal  area  densitv  Xf  =  0  0625 

Canopy  height 

hc 

0.05  m 

Bulk  drag  parameter 

CdAhc 

0.025 

Friction  velocity 

u*  0 

0.0048  m  s'1 

Mean  wind  at  z  =  hc 

U{hc)/u*0 

8.8 

Equilibrium  standard  deviations 

(o-u,  O y  ,  /^*Q 

2.1,  1.6. 

1.2 

TKE  at  z  =  hc 

k(hc)/u*  o 

3.78 

Effective  pressure  gradient 

(■ hc/ul0)OPe/dx 

-0.16 

Frontal 

area  densitv  A f  =  0.16 

Canopy  height 

hc 

0.05  m 

Bulk  drag  parameter 

CdAhc 

0.078 

Friction  velocity 

U* 0 

0.0051  m  s  1 

Mean  wind  at  z  —  hc 

U(hc)/u*  o 

6.2 

Equilibrium  standard  deviations 

,  O'y  ,  CT  w  )/u*q 

2.1,  1.7, 

1.2 

TKE  at  z  —  hc 

k(hc)/ul0 

3  63 

Effective  pressure  gradient 

C hc/ul0)dPe/dx 

-0.22 

Frontal 

area  density  A  f  =  0  44 

Canopy  height 

hc 

0.05  m 

Bulk  drag  parameter 

CdAhc 

0.24 

Friction  velocity 

u*  o 

0.0025  in  s'1 

Mean  wind  at  z  =  hc 

U(hc)/u*  o 

5.14 

Equilibrium  standard  deviations 

(^Ul  C^Uj)/^*0 

2.1,  1.7, 

1.2 

TKE  at  z  =  hc 

4.96 

Effective  pressure  gradient 

{hc/ulQ)0Pe/0x 

-0.17 

Raupach  et  al.  [22]  reported  profiles  of  mean  wind  and  turbulence  statistics  in  the  Tombstone 
canopy  which  was  an  array  of  vertical  bars  (flat  plates)  60  mm  high,  10  mm  wide,  and  1 
mm  thick  arranged  in  a  staggered  (diamond)  pattern  with  60  mm  cross-stream  and  44  mm 
streamwise  spacing.  This  obstacle  arrangement  resulted  in  a  frontal  area  index  of  0  23. 
Velocity  statistics  were  measured  with  a  specially  developed  3-wire  probe  [23]  in  a  wind 
tunnel.  The  parameters  that  characterize  this  canopy  are  summarized  in  Table  6.  Because 
only  two  components  of  the  instantaneous  velocity  (namely,  u  and  w)  were  measured  in  this 
experiment,  the  cross-stream  velocity  standard  deviation  av  and  the  TKE  k  were  obtained 
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assuming  that  =  cruaw  (viz.,  the  cross-stream  standard  deviation  of  velocity  fluctuations 
is  the  geometric  mean  of  the  streamwise  and  vertical  velocity  fluctuations) 


Table  4.  Parameters  of  the  staggered  arrays  of  cubical  obstacle  experiments. 


Property 

Symbol 

Value 

Frontal  area  density  Xp  —  0.0625 
Canopy  height  hc 

0.05  m 

Bulk  drag  parameter 

CdAhc 

0.042 

Friction  velocity 

u+o 

0.0053  m  s"1 

Mean  wind  at  z  =  hc 

U(hc)/u*0 

72 

Equilibrium  standard  deviations 

(o'n  ,  (7  v  .  (Jxu  ) /ti  +  Q 

2.0,  1  6,  1.2 

TKE  at  z  =  hc 

k(hn)/u% 

3.49 

Effective  pressure  gradient 

{hc/u%)dPe/dx 

-0  16 

Frontal 

Canopy  height 

area  density  Xp  =  0  16 

hc 

0.05  m 

Bulk  drag  parameter 

CdAhc 

0 11 

Friction  velocity 

^*0 

0.0062  m  s-1 

Mean  wind  at  z  —  hc 

U(hc)/u, o 

49 

Equilibrium  standard  deviations 

(7V ,  <7w^f  U* o 

2  0,  1  6,  1  2 

TKE  at  z  =  hc 

k(hc)/u*0 

2.94 

Effective  pressure  gradient 

(hc/u'i0)dPe/dx 

-0.24 

Frontal 

Canopy  height 

area  density  Xp  —  0  44 

he 

0.05  in 

Bulk  drag  parameter 

CdAhc 

0.31 

Friction  velocity 

u*0 

0.0026  in  s'1 

Mean  wind  at  z  =  hc 

U(hc)/u*  o 

5.4 

Equilibrium  standard  deviations 

(0*1x1  0"ur)/^*O 

2.2.  1  7,  1.2 

TKE  at  z  =  hc 

k{hc)ful  0 

6.16 

Effective  pressure  gradient 

-0.14 
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Table  5.  Parameters  of  the  square  and  staggered  arrays  of  billboard  obstacle  experiments 
for  frontal  area  density,  Xp  =  0  16. 


Property  Symbol  Value 


Canopy  height 

Bulk  drag  parameter 

Friction  velocity 

Mean  wind  at  z  —  hc 

Equilibrium  standard  deviations 

TKE  at  z  =  hc 

Effective  pressure  gradient 

Canopy  height 

Bulk  drag  parameter 

Friction  velocity 

Mean  wind  at  z  —  hc 

Equilibrium  standard  deviations 

TKE  at  z  —  hc 

Effective  pressure  gradient 


Square  Arrav 

hc 

0  05  m 

CdAhc 

0.35 

0.0082  m  s~ 

C  {he)  /  i^*o 

34 

(u"ui 

2.0,  1  6.  1.2 

k(hc)/ul0 

3.40 

(hc/u%)dPe/dx 

-0.22 

Staggered  Arrav 

hc 

0.05  m 

CdAhc 

0.23 

u*  o 

0.0080  m  s~ 

C  {hc )  j  u*o 

3.2 

{& U  1  &V)  &w  )  i U'* 0 

1.9,  1.5,  1.2 

k{hc)/ul o 

3.13 

(hc/ul0)dPe/dx 

-0.24 

Table  6.  Parameters  of  the  Tombstone  canopy  experiment 


Property 

Symbol 

Value 

Canopy  height 

he 

0.06  m 

Bulk  drag  parameter 

CdAhc 

0.31 

Friction  velocity 

u*  o 

1.1  m  s-1 

Mean  wind  at  z  =  hc 

£/  ( /i(. )  /n  *o 

3.09 

Equilibrium  standard  deviations 

(n u  i  ® v  t 

2.0,  1.5,  1  14 

TKE  at  z  —  hc 

k{hc)/ul o 

3.8 

Effective  pressure  gradient 

(hc/u%)dPe/dx 

-0.23 
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Model  comparisons 


The  numerical  method  used  as  the  basis  for  this  work  is  described  m  detail  m  [10]  Briefly, 
the  governing  fluid  equations  for  both  the  basic  and  alternative  canopy  windflow  models 
were  solved  using  a  finite-volume  discretization  along  with  a  staggered  grid  arrangement 
Over  this  grid  system,  the  governing  fluid  flow  equations  are  discretized  and  this  yields  a 
set  of  equivalent  discrete  algebraic  equations.  These  algebraic  equations  are  then  solved  by 
an  iterative  procedure,  subject  to  the  appropriate  boundary  conditions  which  specify  the 
problem.  The  SIMPLEC  algorithm  [24]  was  used  to  couple  pressure  and  velocity  on  the 
staggered  grid,  and  a  line-by-line  iterative  procedure  was  applied  to  solve  the  discretized 
system  of  algebraic  equations  over  the  computational  domain  Because  all  the  model  simu¬ 
lations  involved  uniform  canopies,  horizontal  homogeneity  is  assumed  leading  to  0/dxx  =  0 
for  i  =  1,2.  All  the  simulations  are  done  in  rotated  coordinates  so  that  the  flow  is  m  the 
xi-direction  (i.e.,  V  =  0). 

Local  equilibrium  lower  boundary  conditions  are  used  for  U  and  k.  In  the  staggered  grid  that 
defines  abutting  control  volumes  for  U  and  /c,  we  have  set  the  lowest  TKE  gridpoints  on  the 
ground.  Consequently,  the  lowest  U  gridpomts  lie  above  the  ground  and  a  “wall  function'’ 
relationship  is  used  to  relate  the  surface  shear-stress  and  the  near-ground  wmdspeed.  In 
particular,  if  the  lowest  U  gridpomts  lie  at  z  —  zp.  then  a  condition  on  the  corresponding 
vertical  momentum  fluxes  to  ground  is  required  In  other  words,  the  boundary  condition  on 
U  (tangential  component  of  U)  is  applied  implicitly  through  the  specification  of  the  shear 
stress  as  u'w'  —  —  u*  where 

u.  =  kUp/Yo.  (— ).  (26) 

V  ^O.s  J 

Here,  k  m  0  4  is  von  Karman's  constant.  Up  =  U(zp).  and  zqs  is  the  substrate  roughness 
length  at  the  ground.  In  practice,  the  value  of  z0s  is  chosen  rather  arbitrarily,  with  simu¬ 
lations  showing  very  little  sensitivity  to  a  wide  range  of  values.  All  the  model  simulations 
presented  in  the  reported  used  z0s  =  0.00001  m.  The  imposition  of  this  lower  boundary  con¬ 
dition  is  consistent  with  a  constant  stress,  logarithmic  mean  profile  at  the  lower  boundary 
(ground),  which  is  the  simplest  consistent  assumption  that  we  can  make  about  conditions  at 
the  ground  surface  below  the  canopy.  Of  course,  strictly  speaking,  these  wall  relationships 
are  not  valid  within  the  canopy  (nor  for  that  matter  in  regions  of  highly-disturbed  flow), 
but  as  the  shear  stress  on  the  ground  beneath  a  canopy  (that  is  not  too  sparse)  is  generally 
small,  its  miscalculation  (by  the  above  relationships)  is  expected  to  carry  negligible  penalty. 
As  the  lower  boundary  condition  on  TKE  in  conjunction  with  Equation  (26),  we  adopted 
the  equilibrium  relationship  kgnc \  =  uljce .  The  alternative  prescription  (dk/dz)z=z0  =  0 
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performed  neither  better  nor  worse. 

We  set  the  uppermost  U  gridpoint  in  a  normal,  full,  control  volume,  bounded  by  the  top 
of  the  computational  domain,  so  that  the  upper-boundary  condition  on  U  was  a  specified 
shear  stress.  Either  we  specified  this  shear  stress  aloft  as  undisturbed  (u*0)  and  constituting 
the  velocity  scale  for  the  simulation  (in  the  case  of  a  field  simulation  such  as  the  Elora  corn 
canopy);  or,  to  properly  account  for  stress  and  TKE  gradients  at  z  >  hc  in  the  incident  flow, 
we  incorporated  an  effective  background  kinematic  pressure  gradient  [cf  Equation  (25)],  and 
specified  the  shear  stress  aloft  as  —u!w\z  =  z?)  =  u*0  +  z?  —  /ic),  where  z t  is  the 
top  of  the  computational  domain  and  dPe/dz  =  (dr/dz)  is  the  measured  stress  gradient 
above  the  canopy  (in  the  case  of  a  wind  tunnel  or  water  channel  simulation).  Since  a  k 
gridpoint  lies  on  the  upper  boundary  of  the  computational  domain,  we  applied  the  Neumann 
boundary  condition  ( dk/Oz)z~ZT  =  0.  This  is  not  quite  correct  for  the  wind  tunnel  or  water 
channel  simulations  because  of  the  TKE  gradients  above  the  canopy,  but  the  alternative 
of  specifying  the  value  of  k  directly  at  the  top  of  the  computational  domain  required  an  a 
priori  knowledge  of  this  value  (which  we  assume  would  not  be  available  generally)  In  all 
our  model  simulations,  we  “turned  off”  the  Coriolis  acceleration  by  setting  /  =  0  because 
the  full  PBL  was  not  required 

Other  rather  familiar  parameters  have  appeared  in  our  canopy  windflow  model  equations,  as 
constants  which  we  want  to  clearly  distinguish  from  the  closure  parameters — to  emphasize 
that  the  success  of  our  simulations  does  not  depend  on  any  flexibility  in  their  specification. 
They  are  d,  and  ce  and  are  determined  by  the  flow  properties  taken  from  within  the  con¬ 
stant  stress  layer  above  the  canopy.  Firstly,  the  turbulence  length  scale  [cf.  Equation  (12)] 
requires  the  assignment  of  Lc 0  which  imposes  a  limiting  value  on  the  magnitude  of  the  outer 
length  scale.  For  a  full  PBL  simulation,  L ^  is  assigned  using  Equation  (13)  However,  for 
numerical  simulations  of  the  wind  tunnel  and  water  channel  data,  we  set  L^/hc  —  1.0  to 
accomodate  the  shallowness  of  the  above-canopy  constant-stress  layer  in  these  flows.  For 
the  field  experiments  (e.g.,  Elora  corn  canopy),  we  assumed  a  deep  constant-stress  layer  for 
the  background  flow  and  set  L ^  ~  oo.  The  numerical  value  of  affects  our  simulated 
profiles  of  [/,  r,  and  k  only  far  above  the  canopy  in  a  region  of  little  interest  to  us.  The  aim 
was  not  to  replicate  the  actual  details  of  the  flow  aloft  (which  for  the  wind  tunnel  and  water 
channel  simulations  would  depend  on  actual  mechanisms  [e.g.,  boundary- layer  trips,  spires, 
upstream  roughness  elements]  used  to  accelerate  the  development  of  the  boundary-layer 
above  the  canopy),  but  merely  to  ensure  that  the  modeling  of  the  disturbed  flow  through 
the  canopy  (which  is  the  region  of  greatest  interest  in  this  study)  was  not  compromised  by 
an  inherently  unrealistic  treatment  aloft.  Secondly,  the  displacement  length  d  required  in 
the  definition  of  the  outer  turbulence  length  scale  A0  [cf.  Equation  (12)]  was  either  assigned 
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the  value  reported  for  the  canopy  (e.g.,  Furry  Hill  and  Elora  corn  canopy  simulations),  or 
estimated  using  a  similarity  relationship  for  d  in  regular  obstacle  arrays  proposed  by  Mac¬ 
donald  et  al.  [25]  (e  g  ,  for  the  wind  tunnel  and  water  channel  ‘'urban*'  canopy  simulations) 
Thirdly,  the  value  of  ce  [cf.  Equation  (11)]  was  determined  for  each  canopy  from  the  mea¬ 
surements  of  the  equilibrium  ratios  cu  —  au/n+ o,  etc.  above  the  canopy  [cf.  Equation  (24)] 
However,  these  equilibrium  ratios  m  flows  above  uniform  canopies  typically  vary  little  from 
one  canopy  to  the  next  Indeed,  typically  ~  2.1, 1.7. 1.2.  implying  ce  ^  0.23. 

This  value  for  ce  could  have  been  used  for  all  the  model  simulations,  with  little  change  in 
the  final  results 


It  only  remained  in  the  model  simulations  to  specify  the  bulk  drag  coefficient  Cd  (note  that 
A  and  hc  are  determined  from  the  geometric  characteristics  of  the  canopy)  Unfortunately, 
there  exists  currently  no  theory  that  allows  one  to  predict  the  bulk  drag  coefficient  of  a 
canopy  from  knowledge  of  the  canopy  geometry  and  the  behavior  of  the  obstacle  elements 
in  isolation  In  our  canopy  wmdflow  model,  we  view  Cd  as  a  data  input.  For  model  evalu¬ 
ation,  the  value  of  Cd  for  a  particular  canopy  was  diagnosed  from  the  measurements  of  the 
mean  wind  and  shear  stress  within  the  canopy  [cf  Equation  (23)]  This  is  in  contrast  to  the 
alternative  of  considering  Cd  as  an  adjustable  model  parameter  whose  value  is  simply  ad¬ 
justed  to  optimise  the  fit  between  modeled  and  measured  mean  profiles.  The  first  approach 
ensures  that  the  Utrue"  rate  of  momentum  extraction  occurs  within  the  canopy  regardless 
of  the  wind  speed 


The  canopy  windflow  model  used  in  this  study  is  a  one-equation  model  that  is  based  on 
an  approximate  transport  equation  for  the  turbulence  kinetic  energy  and  on  the  use  of  the 
Boussinesq  linear  eddy-viscosity  approximation  to  compute  the  Reynolds  stress  tensor  as  the 
product  of  an  eddy  viscosity  and  the  mean  strain-rate  tensor  As  a  further  approximation, 
we  assumed  that  the  noimal  stress  gradients  play  a  much  smaller  role  than  the  shear  stress 
gradients  in  the  determination  of  the  mean  flow,  so  we  neglected  the  normal  stress  gradient 
contributions  to  the  mean  momentum  equations,  replacing  them  instead  with  small  artificial 
diffusion  terms.  These  diffusion  terms  were  added  simply  to  ensure  numerical  stability.  The 
artificial  viscosity  Ka  used  to  characterize  these  diffusion  terms  is  very  small.  Generally, 
for  the  current  model  simulations,  the  artificial  viscosity  was  specified  as  Ka  —  0.01hcu*o 
and  had  insignificant  effect  on  the  simulations  other  than  to  ensure  numerical  stability  For 
model  comparisons,  the  velocity  variances  (normal  stresses)  cr^,  and  cr^  were  derived 
from  the  calculated  TKE,  assuming  equilibrium  partitioning,  i.e.  o\  —  cec^/c,  =  cec^fc. 
and  =  cec^fc.  Note  that  this  reduces  to  the  special  instance  of  isotropic  partitioning  of 
the  TKE  if  and  only  if  cu  =  cv  =  cw ,  which  is  generally  not  the  case  for  canopy  flows. 
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Figures  1,2,  and  3  compare  modeled  and  measured  profiles  of  t/,  r,  and  k  for  the  aeroelastic 
plant  canopy  used  in  the  Furry  Hill  experiment  In  this  and  subsequent  comparisons, 
the  solid  and  long  dashed  curves  correspond  to  predictions  obtained  from  the  basic  and 
alternative  canopy  windflow  models,  respectively  Vertical  axes  are  scaled  by  hc ,  while  mean 
velocity  and  turbulence  profiles  are  scaled  by  u*q  and  u*0,  respectively.  The  basic  windflow 
model  (viz.,  simulations  using  c  =  a*  =  l,/i  =  02,  all  other  inputs  [d:  ce,  Cj ,  A]  provided 
unambiguously  by  the  Furry  Hill  observations)  provides  excellent  agreement  with  the  Furry 
Hill  data,  although  it  must  be  emphasized  that  the  closure  constants  here  were  optimized 
with  respect  to  this  data  set.  The  quantitative  agreement  between  model  and  observations 
for  the  alternative  windflow  model  (viz.,  c  —  p  —  \  with  a  model  [cf.  Equation  (22)]  for  the 
form  drag  sink  in  the  TKE  transport  equation  that  does  not  involve  the  closure  constant 
a*)  is  not  as  good  as  that  for  the  basic  windflow  model,  but  then  again  the  closure  constants 
in  the  alternative  model  have  not  been  optimized  with  respect  to  the  Furry  hill  data. 

The  better  agreement  of  the  basic  model  predictions  and  the  Furry  Hill  observations  relative 
to  the  alternative  model  predictions  arises  primarily  from  the  specification  fi  =  0.2  Recall 
that  the  turbulent  diffusivity  Kk  =  l^K  of  TKE  parameterizes  the  sum  (Tt+Tp)  of  turbulent 
(Tt)  and  pressure  (Tp)  transport  of  TKE.  Also,  recall  that  on  theoretical  grounds,  Lumley 
[26]  suggested  (pf  is  the  kinematic  pressure  fluctuations) 

P'uj  ~  -0-2 u'ufcuj^. 

Direct  numerical  simulations  conducted  by  Rogers  and  Moser  [27]  for  a  plane  mixing  layer 
supported  Lumley’s  suggestion.  The  crucial  point  is  that,  if  \Tt  +  rp|/|Tt|  is  systematically 
different  in  the  canopy  than  in  the  plane  mixing  layer,  then  there  is  no  reason  why  the 
conventional  specification  /i  =  1  (as  recommended  in  [16,17])  should  obtain.  However,  the 
evidence  seems  to  suggest  the  contrary.  Raupach  et  al  [13]  showed  that  the  plane  mixing 
layer  flow  is  a  useful  archetype  for  turbulence  in  the  vicinity  of  a  canopy,  and  used  this 
plane  mixing  layer  analogy  to  provide  a  '‘universal*’  view  of  canopy  turbulence  that  appears 
to  be  consistent  over  a  wide  range  of  canopies.  Furthermore,  the  TKE  balance  with  /i  =  0.2 
results  in  a  smaller  turbulent  transport  term  than  Brunet  et  al  [15]  observed  (viz.,  the 
specification  p  =  0.2,  which  implies  that  a  larger  fraction  of  Tt  than  20%  is  negated  by  Tp 
in  the  canopy  [contrary  to  Lumley’s  suggestion],  is  not  consistent  with  these  measurements 
of  the  turbulent  transport  term).  However,  the  specification  (i  =  0.2  was  kept  in  the  basic 
windflow  model  simply  because  this  specification  resulted  in  a  better  match  to  the  observed 
k(z)  profile  in  the  Furry  Hill  experiment. 

Figures  4,  5,  and  6  display  basic  and  alternative  canopy  model  predictions  of  {/,  r,  and 
fc,  respectively,  for  the  Elora  corn  canopy.  Both  the  basic  and  alternative  models  provide 
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reasonably  good  simulations  of  the  field  measurements  The  alternative  model  provides 
a  slightly  better  agreement  for  the  mean  velocity  profile  than  the  basic  model,  with  the 
converse  being  true  for  the  turbulence  profiles.  For  this  canopy,  the  use  of  a  local  drag 
coefficient  diagnosed  from  the  measured  mean  wind  and  shear  stress  within  the  canopy  as 

resulted  in  a  better  agreement  of  model  and  observations  (for  both  the  basic  and  alternative 
models),  but  because  such  detailed  information  on  the  local  drag  coefficient  C^(z)  would 
usually  be  unavailable,  we  have  not  pursued  this  approach  any  further  and  have  performed 
all  simulations  using  a  bulk  drag  coefficient  CT  With  this  understanding,  it  is  expected 
that  the  model  simulations  will  represent  a  less  perfect  reproduction  of  the  data  than  could 
be  achieved,  but  this  is  ail  accepted  cost  of  keeping  the  system  simple. 

Figures  7  and  8  show  comparisons  of  model  predictions  and  observations  of  the  mean  velocity 
(U)  and  turbulence  (r,  A\  cr2)  for  a  square  array  of  cubical  obstacles  with  a  frontal  area  index 
of  0.0625  Both  the  basic  and  alternative  models  are  seen  to  provide  very  good  predictions 
for  the  mean  velocity  U  and  shear  stress  r  for  this  canopy  However,  it  is  observed  that  the 
basic  model  underestimates  k  and  cr2  deep  in  the  canopy.  The  alternative  model,  on  the 
other  hand,  provides  a  much  better  estimate  of  k  and  cr2  within  the  canopy  Fortunately, 
despite  the  fact  that  the  basic  model  underestimates  k  within  the  canopy,  its  predictions 
for  U  and  r  are  basically  as  good  as  those  of  the  alternative  model.  This  is  because  the 
eddy  viscosity  is  proportional  to  A;1/2,  rendering  U  and  r  remarkably  insensitive  to  wide 
variations  in  the  treatment  of  k. 

The  measured  and  modeled  mean  wind  and  turbulence  profiles  for  the  water  channel  square 
array  of  cubical  obstacles  with  a  frontal  area  index  of  0  16  are  displayed  in  Figures  9,  10. 
and  11.  The  overall  conformance  of  model  and  observations  for  the  basic  and  alternative 
canopy  models  is  similar  to  that  for  the  square  cubical  array  with  A p  —  0.0625  It  is 
apparent  from  Figures  10  and  11  that  the  equilibrium  partitioning  of  TKE  to  derive  the 
three  turbulent  velocity  variances  a2,  cr2,  and  seems  to  work  well  for  this  example  The 
equilibrium  partitioning  of  k  used  here  for  all  vertical  levels  2  followed  the  measured  pattern 
of  the  partioning  of  k  at  the  canopy  top  (which  has  been  summarized  in  Table  3  under  the 
property  of  equilibrium  standard  deviations).  Indeed,  the  magnitude  of  the  three  velocity 
variances  is  seen  to  compare  favorably  with  the  observations.  In  addition,  it  is  noted  that 
because  the  turbulent  velocity  variances  are  constructed  essentially  by  applying  scaling 
factors  to  the  associated  model  profile  for  k  (with  the  scaling  factors  determined  from  the 
values  of  ce,  cu ,  cv,  and  cw  diagnosed  from  the  flow  properties  above  the  canopy),  all  these 
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profiles  have  necessarily  the  same  shape.  Unfortunately,  given  the  scatter  in  the  turbulence 
data  for  this  example,  it  cannot  be  determined  how  accurate  this  shape  assumption  is. 

Next,  we  consider  the  comparison  of  model  results  with  observations  of  U ,  r.  k,  and  a\ 
exhibited  m  Figures  12  and  13  for  the  square  array  of  cubical  obstacles  with  a  frontal 
area  index  of  0  44.  Although  both  the  basic  and  alternative  windflow  models  agree  with 
measurements  of  U  within  the  canopy  (z  <  hc),  its  overall  conformance  with  the  observed 
U  profile  above  the  canopy  is  poor.  In  particular,  both  models  underpredict  severely  the 
magnitude  of  the  mean  wind  shear  at  the  top  of  the  canopy,  and  this  underprediction  leads 
to  the  poor  agreement  of  modeled  and  observed  U  above  the  canopy  This  case  corresponds 
to  the  skimming  flow  regime  where  the  flow  in  the  canopy  is  vertically-decoupled  from  the 
flow  aloft,  and  the  windflow  models  used  here  do  not  seem  to  handle  this  case  well  (at 
least  with  respect  to  the  prediction  of  the  mean  velocity  profiles).  However,  despite  this 
deficiency  in  our  model  predictions,  overall  agreement  of  predictions  with  the  turbulence 
profiles  (e.g.,  shear  stress,  normal  stresses,  and  TKE)  is  very  good.  Interestingly,  the 
alternative  windflow  model  (which  provides  a  relatively  good  prediction  of  the  turbulence 
level  within  the  square  array  of  cubical  obstacles  over  the  range  of  frontal  area  indices 
considered)  predicts  that  the  Reynolds  stresses  (and,  the  TKE  which  is  simply  the  trace  of 
the  Reynolds  stress  tensor)  diminish  from  the  top  of  the  canopy  downwards  exponentially 
for  A p  =  0  44  (cf.  Figures  12  and  13),  whereas  this  rate  of  decrease  is  roughly  linear  for 
the  smaller  frontal  area  densities  (cf.  Figures  7  to  11  for  A^  =  0.0625  and  0  16)  The  data 
appear  to  be  consistent  with  this  model  prediction,  although  measurements  deep  within 
the  square  cube  array  with  XF  =  0.44  are  not  available  for  the  full  confirmation  of  this 
prediction  owing  to  measurement  constraints  imposed  by  the  small  separation  between  the 
obstacles  for  this  case. 

Figures  14,  15,  and  16  compare  the  modeled  (obtained  from  the  basic  and  alternative 
canopy  windflow  models)  and  measured  vertical  profiles  of  t/.  r,  fc,  cr£,  <7^,  and  for  the 
staggered  array  of  cubical  obstacles  with  a  frontal  area  index  of  XF  —  0.0625.  Similarly,  the 
measured  and  modeled  mean  wind  and  turbulence  profiles  for  the  staggered  array  of  cubical 
obstacles  with  a  frontal  area  density  of  0.16  are  shown  in  Figures  17,  18,  and  19.  Finally, 
Figures  20,  21,  and  22  compare  the  mean  wind  and  turbulence  profiles  with  solutions  from 
the  basic  and  alternative  windflow  models  for  the  staggered  array  of  cubical  obstacles  with 
A p  =  0.44.  On  examination  of  these  figures,  it  is  apparent  that  the  overall  conformance  of 
model  predictions  and  observations  for  the  staggered  cube  arrays  is  similar  to  that  for  the 
corresponding  square  cube  arrays. 

Figures  23,  24,  and  25  show  vertical  profiles  of  mean  wind  speed  (U)  and  turbulence  (r,  fc, 
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o\,  &y ,  cn  )  for  the  square  array  of  billboard-shaped  obstacles  with  a  frontal  area  density 
of  0.16  Numerical  simulations  obtained  from  the  basic  and  alternative  canopy  windflow 
models  are  compared  with  the  experimental  data.  Figure  24  shows  reasonable  agreement 
between  the  modeled  and  observed  profiles  for  U  and  r.  However,  unlike  the  cases  for 
the  square  and  staggered  arrays  of  cubical  obstacles,  the  TKE  and  the  turbulent  velocity 
variances  within  the  canopy  are  seen  to  be  underpredicted  by  both  the  basic  and  alternative 
windflow  models.  Unlike  the  square  cube  array  with  Xp  =  0.16,  the  square  billboard  array 
with  the  same  frontal  area  index  possesses  a  much  smaller  plan  area  density.  This  implies 
that  the  square  billboard  array  posseses  much  more  “free-air*’  space  and,  hence,  experiences 
a  greater  exposure  to  the  wind  than  the  square  cube  array.  A  pecularity  of  the  flow  over 
the  billboard  canopy  is  that  the  wind  profile  aloft  is  more  nearly  logarithmic  in  z  than 
in  z  —  d\  viz.,  the  displacement  length  d  %  0,  which  is  consistent  with  predictions  from  a 
similarity  relation  for  d  (e.g.,  Macdonald  et  al.  [25])  In  view  of  the  amount  of  “free- air” 
space  within  the  billboard  canopy,  it  would  seem  that  a  reasonable  procedure  here  would  be 
to  eliminate  the  form-drag  sink  term  contribution  in  the  TKE  transport  equation  (i  e.,  set 
€fd  =  0  and  model  the  viscous  dissipation  as  e  =  max(ccc,  Cfd)  —  ecc)  Figures  24  shows  that 
eliminating  the  form  drag  sink  for  the  TKE  (in  the  alternative  windflow  model)  increases 
k  within  the  billboard  canopy  (Short  dashed  curve)  to  a  level  that  is  more  representative  of 
the  observations,  with  the  concomitant  improvement  m  simulation  of  the  turbulent  velocity 
variances  (cf.  Figures  24  and  25)  Again,  it  should  be  noted  that  these  different  variations 
in  the  treatment  of  k  have  little  effect  on  the  predictions  of  U  and  r. 

The  measured  and  modeled  mean  wind  and  turbulence  profiles  for  the  staggered  array  of 
billboard-shaped  obstacles  with  a  frontal  area  index  of  0.16  are  displayed  in  Figures  26,  27, 
and  28.  Again,  it  is  seen  that  the  elimination  of  the  form  drag  sink  in  the  TKE  transport 
equation  provides  a  better  simulation  of  the  turbulence  levels  (eg.,  k,  cq;,  and  <r^) 
within  a  billboard-shaped  obstacle  canopy.  Given  the  scatter  in  the  measurements  of  the 
turbulence  profiles,  it  is  seen  that  equilibrium  partitioning  of  the  TKE  provides  a  reasonably 
good  reproduction  of  the  three  normal  stresses. 

Finally,  Figures  29,  30,  and  31  show  observed  profiles  of  C7,  t,  and  A;,  respectively,  for  the 
Tombstone  canopy  in  comparison  with  numerical  simulations  using  the  basic  and  alternative 
windflow  models  as  well  as  a  modification  of  the  alternative  windflow  model  with  the  TKE 
form  drag  sink  term  removed.  The  model  profiles  for  U  and  r  are  rather  good.  However, 
unlike  the  preceding  staggered  billboard  canopy,  the  model  prediction  for  k  with  the  form 
drag  sink  in  the  TKE  eliminated  now  overestimates  k  within  the  Tombstone  canopy.  In  this 
case,  the  alternative  windflow  model  (with  the  form  drag  term  retained  in  the  TKE  transport 
equation)  provides  the  best  conformance  with  the  observations  of  k.  This  is  suprising  in  view 
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of  the  preceding  results  for  the  staggered  billboard  canopy  This  is  because  the  Tombstone 
canopy  studied  by  Raupach  et  al.  [22]  is  essentially  a  staggered  array  of  billboard-shaped 
obstacles  with  frontal  area  density  0.23  and,  hence,  is  very  similar  to  the  staggered  billboard 
array  with  Ap  =  0.16  investigated  by  Macdonald  et  al.  [21]  In  particular,  both  these 
obstacle  arrays  correspond  to  near  maximally  rough  flow  (i.e.,  the  roughness  length  z0 
achieves  a  maximum  with  frontal  area  densities  Ap  «  0.15-0.2),  with  the  relationships 
between  velocity  statistics  being  very  similar  the  normalized  mean  wind  speed  at  canopy 
top,  U(hc)/u* o  —  3  1  for  the  Tombstone  canopy  compared  with  U(hc)/u* o  =  3.2  for  the 
billboard  canopy;  and,  the  normalized  TKE  at  canopy  height.  k(hc)/v%0  =  3.8  for  the 
Tombstone  canopy  which  is  close  to  the  value  k{hc)/ul0  =  3  2  for  the  billboard  canopy. 
However,  despite  these  similarities,  it  is  seen  that  the  measured  TKE  within  the  Tombstone 
canopy  is  systematically  smaller  in  value  than  that  observed  in  the  staggered  billboard 
canopy  Furthermore,  the  TKE  diminishes  roughly  linearly  within  the  staggered  billboard 
canopy  (cf  Figure  27),  whereas  the  TKE  in  the  Tombstone  canopy  exhibits  a  nonlinear 
decrease  (1  e.,  exponential-like  decrease)  within  the  canopy  It  is  unknown  whether  these 
seemingly  differing  estimates  of  the  TKE  in  these  two  rather  similar  canopies  is  the  result 
of  measurement  error  or  is  a  physical  effect  arising  from  the  slight  differences  in  the  two 
canopies. 
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Discussion  and  conclusions 


The  basic  and  alternative  canopy  windfiow  models  were  applied  to  simulate  flows  in  a  wide 
range  of  canopies  (both  vegetative  and  urban).  The  models  were  assessed  by  comparison 
with  available  wind  tunnel,  water  channel,  and  field  data,  where  we  have  restricted  our 
attention  to  experiments  providing  measurements  of  both  the  mean  flow  and  the  turbulence 
(shear  stress,  TKE  and/or  normal  stresses)  A  general  and  consistent  procedure  was  used 
for  all  the  model  simulations.  We  use  wall  functions  to  impose  local  equilibrium  boundary 
conditions  at  the  ground,  and  specify  either  a  shear  stress  (for  a  constant  stress  layer  host 
flow)  or  a  pressure  gradient  aloft  (for  a  wind  tunnel  or  water  channel  boundary-layer  host 
flow)  as  appropriate  Our  boundary  conditions  are  otherwise  neutral  in  the  sense  that 
they  do  not  mould  the  profiles  of  £/,  r,  and  k  We  adopt  the  observed  equilibrium  ratios 
cu  —  cru/u* o,  cv  ~  (7v/u+ o,  and  cw  =  <rU)/u¥ o  above  the  canopy  (and,  which  effectively  sets 
the  value  of  the  physical  parameter  ce).  In  any  event,  these  quantities  do  not  vary  greatly 
from  canopy  to  canopy,  or  from  one  equilibrium  layer  to  another  and,  as  such,  one  common 
set  of  values  could  have  been  adopted  for  these  quantities  with  no  appreciable  differences 
in  the  final  results.  We  use  measured  CdA  to  ensure  that  m  the  canopy,  whatever  the 
wind  speed,  the  true  rate  of  momentum  extraction  ensues.  From  these  inputs  and  with  the 
closure  constants  in  eithei  the  basic  or  alternative  windfiow  models  fixed,  we  obtained  quite 
acceptable  predictions  of  profiles  of  mean  wind,  shear  stress,  and  TKE  for  almost  the  entire 
range  of  canopies  investigated  It  was  found  that  the  alternative  windfiow  model  provided 
the  best  overall  conformance  for  a  fairly  wide  range  of  canopies.  The  basic  windfiow  model 
has  the  tendency  to  underpredict  the  TKE  within  the  canopy,  although  this  underprediction 
does  not  appear  to  affect  its  reasonably  good  predictions  of  U  and  r  However,  both  the 
basic  and  alternative  windfiow  models  were  found  to  underpredict  severely  the  magnitude 
of  the  mean  wind  shear  at  the  canopy  top  for  “urban'’  canopies  at  the  highest  frontal  area 
density  [0  44]  investigated  Hence,  it  is  suggested  that  the  alternative  windfiow  model  can 
be  used  probably  unmodified  (viz.,  without  resorting  to  tuning  of  any  model  constants)  for 
the  simulation  of  a  wide  range  of  canopy  flows,  provided  the  frontal  area  densities  of  the 
“urban”  canopies  are  not  so  large  as  to  provoke  skimming  flow.  The  equilibrium  partitioning 
of  TKE  for  the  estimation  of  the  individual  turbulent  velocity  variances  (normal  stresses) 
appears  to  provide  reasonably  good  estimates  for  these  quantities. 

How  critically  have  we  compared  model  results  with  data?  To  address  this  question,  the 
following  points  are  pertinent:  (1)  we  chose  maximal  simplicity  in  the  flow  modeling  to 
minimize  the  introduction  of  flexible  parameters  (closure  constants);  (2)  we  insisted  on 
using  a  constant  bulk  drag  coefficient  Cd  within  the  canopy  that  was  diagnosed  (and,  not 
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adjusted  to  optimize  the  model  fit  to  data)  from  the  within-canopy  mean  wind  and  shear 
stress  profiles,  with  the  implicit  assumption  that  detailed  information  on  a  height-dependent 
Cd  would  usually  be  unavailable;  (3)  free  parameters  (e.g  .  ce,  cu,  d,  etc.)  required  for 
the  model  were  determined  from  the  physical  properties  of  the  flow,  and  we  explicitly 
acknowledged  the  role  of  these  constants  in  adjusting  model  curves  towards  the  data;  (4) 
we  used  a  single  normalizing  velocity  scale  (u*o),  so  as  to  avoid  imparting  any  misleading 
appearance  of  quality  of  the  simulations;  and,  (5)  we  displayed  modeled  and  observed  flow 
quantities  so  as  to  clearly  display  differences.  In  summary,  although  the  success  of  the 
simulations  is  not  automatic,  and  depends  on  the  experience  of  the  user  to  apply  the  model 
appropriately,  making  necessary  judgements  especially  where  important  data  are  uncertain, 
it  remains  that  given  a  set  of  inputs,  the  models  used  here  are  completely  objective  and 
reproducible 

Undoubtedly,  the  most  important  user  input  parameter  to  the  wmdflow  model  is  the  bulk 
drag  coefficient  Cd-  For  all  practical  purposes,  the  most  important  distinction  between 
the  simulations  of  the  various  canopies  lies  in  the  inclusion  of  the  correct  C^.  How  to 
determine  an  appropriate  a  prion  value  of  Cd  for  a  given  arbitrary  canopy  to  use  as  a  data 
input  to  the  model  is  an  outstanding  problem  that  needs  to  be  resolved.  Unfortunately, 
given  our  current  undertanding  of  canopy  flows,  it  is  not  possible  to  predict  the  drag  on 
the  canopy  from  knowledge  of  the  obstacle  geometry  and  the  behavior  of  the  obstacles 
in  isolation.  This  represents  a  major  weakness  in  the  practical  application  of  the  canopy 
windflow  models  proposed  here  Aside  from  the  specification  of  the  data  input  C^,  perhaps 
the  most  uncertain  aspects  of  the  implied  set  of  equations  defining  the  windflow  model  is 
the  treatment  of  TKE  sources  and  sinks  in  the  canopy.  The  basic  and  alternative  wmdflow 
models  used  in  the  present  paper  represent  but  two  different  methods  for  the  treatment 
of  the  TKE  sources  (wake  production)  and  sinks  (form  drag  dissipation).  Our  practice  of 
writing  c  =  max(ecc,efd)  is  inelegant,  and  possibly  not  even  advantageous.  Certainly,  there 
is  a  fundamental  uncertainty  about  how  best  to  treat  these  dissipation  mechanisms  How 
best  to  partition  (i.e.,  spectrally  divide)  the  available  TKE  (from  the  various  sources  such 
as  shear  production,  wake  production,  turbulent  transport)  into  free-air  dissipation  and 
the  extra  canopy  dissipation  arising  from  its  interaction  with  the  obstacles  comprising  the 
canopy  remains  to  be  resolved. 

In  summary,  the  turbulence  closures  used  in  the  windflow  models  here,  which  is  perhaps 
something  like  the  minimum  level  of  closure  required  to  calculate  flows  reasonably  well 
within  a  wide  range  of  canopies,  captures  quite  well  features  of  the  mean  wind  and  tur¬ 
bulence  in  a  canopy.  This  is  encouraging  because  without  any  alteration  of  the  closure 
constants  for  the  alternative  windflow  model  (and  which  have  not  actually  been  calibrated 
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against  any  particular  experimental  data  set  for  canopy  flows),  we  have  achieved  reason¬ 
ably  good  simulations  of  laboratory  and  field  measurements.  This  is  perhaps  surprising,  for 
the  closure  makes  little  more  than  a  token  recognition  of  the  complexity  of  canopy  flows. 
Indeed,  its  sole  concession  to  recent  thought  on  the  mechanisms  of  canopy  flow  is  to  invoke 
the  plane  mixing  layer  analogy  to  provide  a  canopy  length  scale  This  seeming  moderate 
success  in  accounting  for  canopy  flows  using  rather  simple  models  proposed  here  may  result 
from  the  fact  that  the  striking  spatial  variations  of  mean  flow  properties  in  a  canopy  can 
largely  be  explained  by  mean  advection  and  turbulent  diffusion  (parameterized  here  using 
a  gradient  diffusion  hypothesis)  in  the  presence  of  average  sources  (rates  of  momentum 
extraction,  kinetic  energy  transformations,  etc  ). 
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Figure  1.  Mean  wind  speed  observed  in  and  above  a  model  aeroelastic  plant 
canopy  in  a  wind-tunnel  boundary-layer  (Furry  Hill  experiment)  compared  with 
the  outcome  of  model  simulations  using  the  basic  (solid  curve)  and  alternative 
(long  dashed  curve)  canopy  windflow  models. 
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Figure  2.  Shear  stress  observed  in  and  above  a  model  aeroelastic  plant  canopy 
in  a  wind-tunnel  boundary-layer  (Furry  Hill  expen'ment)  compared  with  the 
outcome  of  mode!  simulations  using  the  basic  (solid  curve)  and  alternative  (long 
dashed  curve)  canopy  windflow  models. 
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Figure  3.  Turbulence  kinetic  energy  observed  in  and  above  a  model  aeroelastic 
plant  canopy  in  a  wind-tunnel  boundary-layer  (Furry  Hill  experiment)  compared 
with  the  outcome  of  model  simulations  using  the  basic  (solid  curve)  and 
alternative  (long  dashed  curve)  canopy  windflow  models. 
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Figure  4.  Mean  wind  speed  observed  in  a  canopy  of  corn  (Elora,  Ontario) 
compared  with  the  outcome  of  model  simulations  using  the  basic  (solid  curve) 
and  alternative  (long  dashed  curve)  canopy  windflow  models. 
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Figure  5.  Shear  stress  observed  in  a  canopy  of  corn  (Elora,  Ontario)  compared 
with  the  outcome  of  model  simulations  using  the  basic  (solid  curve)  and 
alternative  (long  dashed  curve)  canopy  windflow  models. 
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Figure  6.  Turbulence  kinetic  energy  observed  in  a  canopy  of  corn  (Elora, 
Ontario)  compared  with  the  outcome  of  model  simulations  using  the  basic  (solid 
curve)  and  alternative  (long  dashed  curve)  canopy  wind  flow  models. 
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Figure  7.  Mean  wind  speed  (left  panel)  and  shear  stress  (right  panel)  observed  in 
and  above  a  square  array  of  cubical  obstacles  with  a  frontal  area  index  of  0.0625 
compared  with  the  outcome  of  model  simulations  using  the  basic  (solid  curve) 
and  alternative  (long  dashed  curve)  canopy  windflow  models. 


DRESTR  2001-058 


39 


P516886.PDF  [Page:  52  of  78] 


Figure  8.  Turbulence  kinetic  energy  (left  panel)  and  streamwise  velocity  variance 
(right  panel)  observed  in  and  above  a  square  array  of  cubical  obstacles  with  a 
frontal  area  index  of  0.0625  compared  with  the  outcome  of  model  simulations 
using  the  basic  (solid  curve)  and  alternative  (long  dashed  curve)  canopy  windflow 
models. 
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Figure  9.  Mean  wind  speed  (left  panel)  and  shear  stress  (right  panel)  observed  in 
and  above  a  square  array  of  cubical  obstacles  with  a  frontal  area  index  of  0.16 
compared  with  the  outcome  of  model  simulations  using  the  basic  (solid  curve) 
and  alternative  (long  dashed  curve)  canopy  windflow  models. 


DRESTR  2001-058 


41 


P516886.PDF  [Page:  54  of  78] 


■st 

"?3 


o 


Figure  10.  Turbulence  kinetic  energy  (left  panel)  and  streamwise  velocity  variance 
(right  panel)  observed  in  and  above  a  square  array  of  cubical  obstacles  with  a 
frontal  area  index  of  0.16  compared  with  the  outcome  of  model  simulations  using 
the  basic  (solid  curve)  and  alternative  (long  dashed  curve)  canopy  windflow 
models. 
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Figure  11.  Cross-stream  velocity  variance  (left  panel)  and  vertical  velocity 
variance  (right  panel)  observed  in  and  above  a  square  array  of  cubical  obstacles 
with  a  frontal  area  index  of  0.16  compared  with  the  outcome  of  model  simulations 
using  the  basic  (solid  curve)  and  alternative  (long  dashed  curve)  canopy  windflow 
models. 


DRES  TR  2001-058 


43 


P516886.PDF  [Page:  56  of  78] 


r/l4 


Figure  12.  Mean  wind  speed  (left  panel)  and  shear  stress  (right  panel)  observed 
in  and  above  a  square  array  of  cubical  obstacles  with  a  frontal  area  index  of  0.44 
compared  with  the  outcome  of  model  simulations  using  the  basic  (solid  curve) 
and  alternative  (long  dashed  curve)  canopy  windflow  models. 
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Figure  13.  Turbulence  kinetic  energy  (left  panel)  and  streamwise  velocity  variance 
(right  panel)  observed  in  and  above  a  square  array  of  cubical  obstacles  with  a 
frontal  area  index  of  0.44  compared  with  the  outcome  of  model  simulations  using 
the  basic  (solid  curve)  and  alternative  (long  dashed  curve)  canopy  windflow 
models. 
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Figure  14.  Mean  wind  speed  (left  panel)  and  shear  stress  (right  panel)  observed 
in  and  above  a  staggered  array  of  cubical  obstacles  with  a  frontal  area  index  of 
0.0625  compared  with  the  outcome  of  model  simulations  using  the  basic  (solid 
curve)  and  alternative  (long  dashed  curve)  canopy  windflow  models. 
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Figure  15.  Turbulence  kinetic  energy  (left  panel)  and  streamwise  velocity  variance 
(right  panel)  observed  in  and  above  a  staggered  array  of  cubical  obstacles  with  a 
frontal  area  index  of  0.0625  compared  with  the  outcome  of  model  simulations 
using  the  basic  (solid  curve)  and  alternative  (long  dashed  curve)  canopy  windflow 
models. 
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Figure  16.  Cross-stream  velocity  variance  (left  panel)  and  vertical  velocity 
variance  (right  panel)  observed  in  and  above  a  staggered  array  of  cubical 
obstacles  with  a  frontal  area  index  of  0.0625  compared  with  the  outcome  of 
model  simulations  using  the  basic  (solid  curve)  and  alternative  (long  dashed 
curve)  canopy  windflow  models. 
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Figure  17.  Mean  wind  speed  (left  panel)  and  shear  stress  (right  panel)  observed 
in  and  above  a  staggered  array  of  cubical  obstacles  with  a  frontal  area  index  of 
0.16  compared  with  the  outcome  of  model  simulations  using  the  basic  (solid 
curve)  and  alternative  (long  dashed  curve)  canopy  windflow  models. 
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Figure  18.  Turbulence  kinetic  energy  (left  panel)  and  streamwise  velocity  variance 
(right  panel)  observed  in  and  above  a  staggered  array  of  cubical  obstacles  with  a 
frontal  area  index  of  0.16  compared  with  the  outcome  of  model  simulations  using 
the  basic  (solid  curve)  and  alternative  (long  dashed  curve)  canopy  windflow 
models. 
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Figure  19.  Cross-stream  velocity  variance  (left  panel)  and  vertical  velocity 
variance  (right  panel)  observed  in  and  above  a  staggered  array  of  cubical 
obstacles  with  a  frontal  area  index  of  0.0625  compared  with  the  outcome  of 
model  simulations  using  the  basic  (solid  curve)  and  alternative  (long  dashed 
curve)  canopy  windflow  models. 
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Figure  20.  Mean  wind  speed  (left  panel)  and  shear  stress  (right  panel)  observed 
in  and  above  a  staggered  array  of  cubical  obstacles  with  a  frontal  area  index  of 
0.44  compared  with  the  outcome  of  model  simulations  using  the  basic  (solid 
curve)  and  alternative  (long  dashed  curve)  canopy  wind  flow  models. 


52 


DRESTR  2001-058 


P516886.PDF  [Page:  65  of  78] 


Figure  21.  Turbulence  kinetic  energy  (left  panel)  and  streamwise  velocity  variance 
(right  panel)  observed  in  and  above  a  staggered  array  of  cubical  obstacles  with  a 
frontal  area  index  of  0.44  compared  with  the  outcome  of  model  simulations  using 
the  basic  (solid  curve)  and  alternative  (long  dashed  curve)  canopy  wind  flow 
models. 
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Figure  22.  Cross-stream  velocity  variance  (left  panel)  and  vertical  velocity 
variance  (right  panel)  observed  in  and  above  a  staggered  array  of  cubical 
obstacles  with  a  frontal  area  index  of  0.44  compared  with  the  outcome  of  model 
simulations  using  the  basic  (solid  curve)  and  alternative  (long  dashed  curve) 
canopy  windflow  models. 
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Figure  23.  Mean  wind  speed  (left  panel)  and  shear  stress  (right  panel)  observed 
in  and  above  a  square  array  of  billboard-shaped  obstacles  with  a  frontal  area 
index  of  0.16  compared  with  the  outcome  of  model  simulations  using  the  basic 
(solid  curve)  and  alternative  (long  dashed  curve)  canopy  windflow  models.  The 
short  dashed  curve  corresponds  to  a  model  simulation  using  the  alternative 
windflow  model  with  the  form  drag  sink  for  TKE  eliminated. 
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Figure  24.  Turbulence  kinetic  energy  (left  panel)  and  streamwise  velocity  variance 
(right  panel)  observed  in  and  above  a  square  array  of  billboard-shaped  obstacles 
with  a  frontal  area  index  of  0.16  compared  with  the  outcome  of  model  simulations 
using  the  basic  (solid  curve)  and  alternative  (long  dashed  curve)  canopy  windflow 
models.  The  short  dashed  curve  corresponds  to  a  model  simulation  using  the 
alternative  windflow  model  with  the  form  drag  sink  for  TKE  eliminated. 
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Figure  25.  Cross-stream  velocity  variance  (left  panel)  and  vertical  velocity 
variance  (right  panel)  observed  in  and  above  a  square  array  of  billboard-shaped 
obstacles  with  a  frontal  area  index  of  0.16  compared  with  the  outcome  of  model 
simulations  using  the  basic  (solid  curve)  and  alternative  (long  dashed  curve) 
canopy  windflow  models.  The  short  dashed  curve  corresponds  to  a  model 
simulation  using  the  alternative  windflow  model  with  the  form  drag  sink  for  TKE 
eliminated. 
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Figure  26.  Mean  wind  speed  (left  panel)  and  shear  stress  (right  panel)  observed 
in  and  above  a  staggered  array  of  billboard-shaped  obstacles  with  a  frontal  area 
index  of  0.16  compared  with  the  outcome  of  model  simulations  using  the  basic 
(solid  curve)  and  alternative  (long  dashed  curve)  canopy  windflow  models.  The 
short  dashed  curve  corresponds  to  a  model  simulation  using  the  alternative 
windflow  model  with  the  form  drag  sink  for  TKE  eliminated. 
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Figure  27.  Turbulence  kinetic  energy  (left  panel)  and  streamwise  velocity  variance 
(right  panel)  observed  in  and  above  a  staggered  array  of  billboard-shaped 
obstacles  with  a  frontal  area  index  of  0.16  compared  with  the  outcome  of  model 
simulations  using  the  basic  (solid  curve)  and  alternative  (long  dashed  curve) 
canopy  windflow  models.  The  short  dashed  curve  corresponds  to  a  model 
simulation  using  the  alternative  windflow  model  with  the  form  drag  sink  for  TKE 
eliminated. 
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Figure  28.  Cross-stream  velocity  variance  (left  panel)  and  vertical  velocity 
variance  (right  panel)  observed  in  and  above  a  staggered  array  of  billboard¬ 
shaped  obstacles  with  a  frontal  area  index  of  0.16  compared  with  the  outcome  of 
model  simulations  using  the  basic  (solid  curve)  and  alternative  (long  dashed 
curve)  canopy  windflow  models.  The  short  dashed  curve  corresponds  to  a  model 
simulation  using  the  alternative  windflow  model  with  the  form  drag  sink  for  TKE 
eliminated. 
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Figure  29.  Mean  wind  speed  observed  in  and  above  a  model  canopy  of 
staggered  flat  plates  with  frontal  area  index  of  0.23  (Tombstone  canopy) 
compared  with  the  outcome  of  model  simulations  using  the  basic  (solid  curve) 
and  alternative  (long  dashed  curve)  canopy  windflow  models.  The  short  dashed 
curve  corresponds  to  a  model  simulation  with  the  alternative  windflow  model 
with  the  form  drag  sink  for  TKE  eliminated. 
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Figure  30.  Shear  stress  observed  in  and  above  a  model  canopy  of  staggered  flat 
plates  with  frontal  area  index  of  0.23  (Tombstone  canopy)  compared  with  the 
outcome  of  model  simulations  using  the  basic  (solid  curve)  and  alternative  (long 
dashed  curve)  canopy  windflow  models.  The  short  dashed  curve  corresponds  to 
a  model  simulation  with  the  alternative  windflow  model  with  the  form  drag  sink 
for  TKE  eliminated. 
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Figure  31.  Turbulence  kinetic  energy  observed  in  and  above  a  model  canopy  of 
staggered  flat  plates  with  frontal  area  index  of  0.23  (Tombstone  canopy) 
compared  with  the  outcome  of  model  simulations  using  the  basic  (solid  curve) 
and  alternative  (long  dashed  curve)  canopy  windflow  models.  The  short  dashed 
curve  corresponds  to  a  model  simulation  with  the  alternative  windflow  model 
with  the  form  drag  sink  for  TKE  eliminated. 
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